Research Article 


Journal of Orthoptera Research 2023, 32(2): 189-200 


Geographic variation in phenotypic divergence between two hybridizing 


field cricket species 


Amy R. Byerty!, CLARA JENCK!, ALEXANDER R. B. Goetz!, David B. WEIssmAN2, Davip A. Gray, CHARLES L. Ross’, 


LUANA S. MAROJA®, ERICA L. LARSON! 


1 Department of Biological Sciences, University of Denver, Denver, CO 80208, USA. 
2 Department of Entomology, California Academy of Sciences, Golden Gate Park, San Francisco, CA 94118, USA. 
3 Department of Biology, California State University, Northridge, CA 91330, USA. 


4 School of Natural Science, Hampshire College, Amherst, MA 01002, USA. 
5 Department of Biology, Williams College, Williamstown, MA 01267, USA. 


Corresponding authors: Amy R. Byerly (amy.byerly@du.edu); Erica L. Larson (erica.larson@du.edu) 


Academic editor: Hojun Song | Received 22 July 2022 | Accepted 23 March 2023 | Published 25 September 2023 


https://zoobank. org/D3283705-EA19-4CB8-9ED6-E113287B2149 


Citation: Byerly AR, Jenck C, Goetz ARB, Weissman DB, Gray DA, Ross CL, Maroja LS, Larson EL (2023) Geographic variation in phenotypic divergence 
between two hybridizing field cricket species. Journal of Orthoptera Research 32(2): 189-200. https://doi.org/10.3897/jor.32.90713 


Abstract 


Patterns of morphological divergence across species’ ranges can pro- 
vide insight into local adaptation and speciation. In this study, we com- 
pared phenotypic divergence among 4,221 crickets from 337 populations 
of two closely related species of field cricket, Gryllus firmus and G. pennsyl- 
vanicus, and their hybrids. We found that these species differ across their 
geographic range in key morphological traits, such as body size and ovi- 
positor length, and we directly compared phenotype with genotype for 
a subset of crickets to demonstrate nuclear genetic introgression, pheno- 
typic intermediacy of hybrids, and essentially unidirectional mitochon- 
drial introgression. We discuss how these morphological traits relate to 
life history differences between the two species. Our comparisons across 
geographic areas support prior research suggesting that cryptic variation 
within G. firmus may represent different species. Our study highlights how 
variable morphology can be across wide-ranging species and the impor- 
tance of studying reproductive barriers in more than one or two transects 
of a hybrid zone. 
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Introduction 


Phenotypic divergence can provide insight into evolutionary 
processes acting across different scales of biological organiza- 
tion. Within a single species, phenotypic divergence can reflect 
differences between environments, population histories, or a 
combination of these factors (Gavrilets et al. 2001, Uyeda et al. 
2009, Runemark et al. 2010, Oneal and Knowles 2013, Jenck et al. 
2020). Phenotypic divergence can signal the possible early stages 


of species differentiation (Wolf et al. 2008, Gonzalez et al. 2011, 
Skoglund et al. 2015) and, in closely related species, can shed light 
on local adaptation and patterns of increasing divergence (Britch 
and Cain 2001, Shaw and Mullen 2011). Most studies of species 
divergence have limited replication across the ranges of a species 
pair, and the specific traits that maintain reproductive barriers be- 
tween species are not always clear (Harrison and Larson 2016). 
Geographically comprehensive surveys of phenotypic divergence 
are much harder (Jiménez and Ornelas 2015, Wang et al. 2017, 
Polly and Wojcik 2019, Moran et al. 2020) but critical if we are 
to understand the origin and maintenance of species’ boundaries. 

The relationship between divergent phenotypic characteristics 
and reproductive barriers is most easily studied in places where 
the ranges of closely related species overlap and heterospecific in- 
dividuals mate and produce offspring (Barton and Hewitt 1985, 
Harrison 1990). In the resulting hybrid zone, as the different spe- 
cies co-exist, compete, and interbreed, phenotypic characteristics 
may be more variable among individuals when compared to the 
pure allopatric populations that lie outside the hybrid zone (Hol- 
lander et al. 2018, Sottas et al. 2018). By comparing the phenotypic 
variation between both conspecific allopatric and sympatric popu- 
lations and between heterospecific populations, it becomes possi- 
ble to examine the potential causes of phenotypic evolution, speci- 
ation, and how those mechanisms lead to the reproductive barriers 
that maintain species boundaries (Shaw and Mullen 2011). 

In this study, we examined the phenotypic divergence be- 
tween two closely related and geographically widespread species 
of North American field crickets, Gryllus pennsylvanicus Burmeister 
1838 and G. firmus Scudder 1902, whose common ancestry dates 
to roughly 200,000 years ago (Willett et al. 1997, Maroja et al. 
2009a). The more northern, inland species, G. pennsylvanicus, is 


Copyright Amy R. Byerly et al. This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 4.0), which permits unrestricted use, 
distribution, and reproduction in any medium, provided the original author and source are credited. 


JOURNAL OF ORTHOPTERA RESEARCH 2023, 32(2) 


190 


broadly distributed throughout the United States, while the more 
southern, coastal species, G. firmus, is restricted to the east coast 
and west into Texas (Alexander 1968, Harrison and Arnold 1982, 
Weissman and Gray 2019). These species form a hybrid zone 
along the eastern front of the Appalachian Mountains (Harrison 
and Arnold 1982), and where they co-occur, they are isolated by 
multiple reproductive barriers. The most striking barrier is a one- 
way incompatibility: G. firmus females mated to G. pennsylvanicus 
males lay few eggs that do not hatch (Harrison 1983, Maroja et al. 
2009b, Larson et al. 2012). These two species are also isolated by 
habitat: G. firmus is often found in sandy habitats and has lighter 
coloration and longer ovipositors that can presumably lay eggs 
deeper in sandy soils (Harrison 1986, Ross and Harrison 2006). 
Also, Gryllus firmus is a larger cricket, though size may vary with 
the length of the growing season (Masaki 1961). In some parts of 
the hybrid zone, G. firmus develops faster and emerges earlier in 
the season, leading to temporal isolation (Harrison 1985). 

These morphological differences have been well characterized 
in a handful of locations within the hybrid zone (e.g., Connecti- 
cut), but whether these morphological traits are consistently dif- 
ferent between G. firmus and G. pennsylvanicus remains an open 
question (Weissman and Gray 2019). When species differences are 
studied in only a few locations, it may be impossible to distin- 
guish species-specific traits from within-species local adaptation. 
Morphological traits, such as lighter color and longer ovipositors, 
may have evolved in specific areas due to habitat selection. Like- 
wise, body size may vary with climate and latitude. This paper pre- 
sents the first geographically comprehensive comparison of G. fir- 
mus and G. pennsylvanicus by combining published and unpub- 
lished morphological datasets for these two species across their 
geographic ranges. Our dataset includes 4,221 crickets from 337 
populations, spanning collections over four decades. We had three 
objectives. First, we quantified morphological divergence within 
and between species across their geographic ranges. Second, for 
populations near the hybrid zone, we tested whether traits that 
distinguish species correlate with ancestry. Finally, we examined 
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the correlation between morphological traits and environmental 
variables across the ranges of these species. In doing so, we aimed 
to gain a greater understanding of how population variation and 
local adaptation contribute to divergence and speciation. 


Materials and methods 


Cricket collections. —We compiled a dataset of 4,221 crickets, the 
majority being G. pennsylvanicus but also G. firmus and their hy- 
brids, from 337 collecting localities (Fig. 1). Crickets were sampled 
throughout the United States and Canada, with the largest collec- 
tions coming from the northeastern United States and the hybrid 
zone. Sampling spanned 40 years (1983-2022), with collections 
performed by A.R. Byerly, E.L. Larson, L.S. Maroja, C.L. Ross, and 
R.G. Harrison. In addition to these previously unpublished mor- 
phological data, we included data from Ross and Harrison (2002), 
Larson et al. (2013), and Weissman and Gray (2019), with the lat- 
ter being the most geographically widespread dataset. We also in- 
cluded morphological data from a newly described cricket species, 
G. thinos Weissman and Gray 2019, , which is closely related to 
G. pennsylvanicus and G. firmus (Gray et al. 2020). We included 
G. thinos to enable us to compare morphological variation within 
G. firmus to that of a closely related species that occupies the same 
habitat but is classified as a separate species. 

We categorized each collecting location as allopatric or sympat- 
ric based on past sampling of the field cricket hybrid zone (Harrison 
and Arnold 1982, Willett et al. 1997, Maroja et al. 2009a, Larson 
et al. 2013a, 2014). Populations in and near the hybrid zone often 
have individuals that are pure G. firmus or pure G. pennsylvanicus, 
but they also have many backcrosses and recent generation hybrids 
(Harrison and Bogdanowicz 1997, Maroja et al. 2009a, Larson et 
al. 2013a, 2014). Because of this, we considered any collecting lo- 
cations that were near the hybrid zone to be “sympatric”. We also 
assigned each collecting location to a geographic region (labeled in 
Fig. 1). These regions, identified using climatological data (Karl and 
Koss 1984), were as follows: central (CTR: IL, IN, KY, MO, OH, TN, 


Fig. 1. Map of North American cricket collecting locations. Allopatric populations of Gryllus firmus are in yellow, G. pennsylvanicus are 
in teal, G. thinos populations are in purple, and sympatric G. firmus and G. pennsylvanicus populations are in red. The size of the circle 
corresponds to the sample size for each location. A. Entire range of collection locations in the United States and Canada; B. Enlarged 
area of densely sampled locations in northeast, central, and southeast United States. 
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WV); east north central (ENC: IA, MI, MN, WI); northeast (NE: CT, 
DE, ME, MD, MA, NH, NJ, NY, PA, RI, VT); northwest (NW: ID, OR, 
WA); south (SO: AR, KS, LA, MS, OK, TX); southeast (SE: AL, FL, 
GA, NC, SC, VA); southwest (SW: AZ, CO, NM, UT); west (WE: CA, 
NV); and west north central (WNC: MT, NE, ND, SD, WY). 

In all cases, crickets were collected by hand and maintained in 
plastic containers with food (cat and rabbit food), water vials, and 
shelter prior to freezing. Most samples were collected as adults, 
but in some cases, crickets were collected as late instar nymphs. 
Nymphs were allowed to mature to the adult stage in the laborato- 
ty before freezing. Most collections were done in August-Septem- 
ber, but some crickets were collected in late July or early October. 


Morphological measurements.—We focused only on traits that were 
measured using the same methods across different studies. Crickets 
were measured for body size, as gauged by either body length, femur 
length, and/or pronotum width. Body length was measured from 
the vertical surface of the face to the tip of the abdomen, straight- 
ening the body when necessary. Pronotum width was measured at 
the widest part of the pronotum. Femur length was measured from 
the proximal to distal end of the hind femur. Female ovipositor 
length was measured from the point of attachment on the abdo- 
men to the distal end of the ovipositor. Because ovipositor length 
varies isometrically with body size (Suppl. material 1: fig. S1), we 
also calculated relative ovipositor length as the length of the ovi- 
positor divided by pronotum width or femur length, depending on 
sample availability. We obtained all measurements using Vernier 
calipers and recorded values to the nearest 0.1 mm. 

For a subset of samples where tegmina were available (31 al- 
lopatric crickets and 437 sympatric crickets), we measured their 
color using a USB4000 spectrophotometer with an Ocean Optics 
PX-2 pulsed xenon lamp and SpectraSuite v2.0 software. We mount- 
ed a probe on a metal stand at a 90° angle 0.7 mm from the surface 
of the tegmina. For each male, we recorded and averaged spectral re- 
flectance for three points near the center of the tegmina. We recorded 
spectral measurements as the percentage of reflected light relative to 
a Spectralon white standard, restricted our analyses to wavelengths 
of 300 700 nm, and used a segmental classification method to es- 
timate brightness, chroma, and hue using CLR v1.1 (Montgomery 
2008). We calculated total brightness (B) as R300 700, which is the 
summed reflectance from 300 nm to 700 nm. We also divided our 
reflectance data into four bins of 100 nm each, calculated the total 
brightness for each bin (Br=600-700, By=500-600, Bg=400-500, 
and Bb=300-400), and then calculated chroma: V(BrBg)?+(ByBb)? 
and hue: arctan|(ByBb)/B]/[(BrBg)/B]. 


Molecular markers.—A subset of the crickets in our dataset was pre- 
viously genotyped for mitochondrial DNA haplotype (N = 1,132, 
Harrison et al. 1987, Harrison and Bogdanowicz 1997, Willett 
et al. 1997, Maroja et al. 2009a, Larson et al. 2013b) and/or 110 
Single Nucleotide Polymorphisms (SNPs) from nuclear genes 
with elevated divergence between G. pennsylvanicus and G. fir- 
mus (N = 559, Larson et al. 2013a, 2014). Mitochondrial DNA 
haplotype was determined by sequencing cytochrome c oxidase 
I, the adjacent tRNA-Leu, and a portion of cytochrome c oxidase 
II (Harrison et al. 1987, Willett et al. 1997). SNPs were identi- 
fied from transcriptomes of male accessory glands from two focal 
populations (Ithaca, NY and Guilford, CT; Andrés et al. 2013) were 
genotyped using the Sequenom MassARRAY platform (Larson et 
al. 2013a, 2014). We used these genotype data to recalculate the 
hybrid index while accounting for hemizygosity for male X-linked 
markers using the methods from Shastry et al. (2021). This was 
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especially important because nearly half of these 110 SNPs are 
located on the X chromosome (Maroja et al. 2015, Gainey et al. 
2018). We defined the hybrid index as the proportion of alleles 
that were inherited from G. firmus (hybrid index = 1; Guildford, CT 
(GUI); Tom’s River, NJ (TOM); and Parksley, MD (MET, a.k.a. PAR 
in Larson et al. 2013a, 2014) and G. pennsylvanicus (hybrid index = 
0; Ithaca, NY (ITH); Scranton, PA (SCR); State College, PA (SCO)). 


Analysis of morphological traits and molecular markers.—All analyses 
were conducted in R v4.1.2 (R Core Team 2020). To manipulate 
the data, we used the R packages dplyr v1.0.6 and tidyverse v1.3.1. 
To plot our data, we used the R packages ggplot2 v3.3.5 and ggpubr 
v0.4.0, and to make our maps, we used Maps v3.3.0. For statistical 
analyses, we used commands from the R packages MASS v7.3-54 
and car v3.0-12. We used the R packages corrplot v0.92 and Hmisc 
v4.5-0 to determine environmental variable correlation. We used 
the R packages AlCcmodavg v2.3-1 and MnMIn 1.43.1 to rank 
models based on Akaike Information Criterion and test models. 

To test for differences in morphological traits between species 
and regions, we used the Kruskal-Wallis test followed by a pair- 
wise Wilcoxon rank sum test (PWRST) to determine differences 
between multiple groups. We chose these non-parametric tests 
because our dataset failed Levene’s test for homogeneity of vari- 
ance. We quantified how well morphological traits could classify 
crickets using a linear discriminant analysis (LDA) on allopatric 
crickets. For all analyses, we present the unadjusted p-values and 
indicate in bold the values that were significant following FDR 
correction (Benjamini and Hochberg 1995). 


Environmental predictors of species distributions.—We tested the re- 
lationships between phenotype and environmental variables that 
we predicted would be important in determining species range or 
local adaptation on two scales: 1) across species ranges and 2) at 
an intermediate scale in a well-characterized region of the hybrid 
zone (Connecticut). Across the species ranges, we used only al- 
lopatric crickets that were most clearly differentiated by morphol- 
ogy, and at the intermediate scale, we used both allopatric and 
sympatric crickets. We focused on the two phenotypes that best 
distinguished the two species and that were quantified in most of 
our samples: ovipositor length and pronotum width. 

We identified 10 environmental variables that might be good 
predictors of species’ distributions based on the natural history 
of these species and prior studies of the field cricket hybrid zone 
(longitude, latitude, elevation, precipitation, minimum tempera- 
ture, maximum temperature, human footprint, and three soil 
characteristics; see Larson et al. 2013b). Elevation, precipitation, 
and temperature data were collected from the PRISM Climate 
Group website (https://prism.oregonstate.edu/). Elevation was 
calculated using an 800-m digital elevation model of the conti- 
nental United States. For each site, we collected precipitation vari- 
ables and minimum and maximum temperatures for the year in 
which each cricket was collected. PRISM data were not available 
for sites in Canada. Soil data were collected from the USDA STATS- 
GO2 soil survey (US sites: https://www.nrcs.usda.gov/wps/portal/ 
nrcs/detail/soils/survey/geo/?cid=nrcs142p2_053629) and_ the 
Soil Landscapes of Canada database (Canada sites: https://sis.agr. 
gc.ca/cansis/nsdb/slc/v3.2/index.html). For a subset of sites in the 
northeastern United States, we used soil data from ISRIC SoilGrids 
(Poggio et al. 2021) due to the smaller spatial scale. These data 
were accessed and compiled using the R package soiIDB v2.6.14. 
We used the following variables: average percent sand, average 
percent clay, and average percent organic matter. Due to the high 
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intercorrelation of soil variables confirmed through correlation 
matrix, we excluded average soil percent silt from further analyses. 
We also obtained spatial data from the Last of the Wild Global 
Human Footprint dataset (version 3), consisting of anthropogenic 
impact measured by population density, land use, and transporta- 
tion access at a 1-km resolution (Venter et al. 2016, 2018). 

We used model selection tests that included these 10 environ- 
mental variables to find the combination of variables that best 
explains morphological variation. We ranked competing models 
using Akaike Information Criterion (AIC), and we reported the 
models with the highest goodness-of-fit. 


Data accessibility —All morphological data and collection site in- 
formation, including GPS coordinates and environmental data 
and scripts, are published in Dryad (doi:10.5061/dryad.jwstqjqdx). 


Results 


Estimates of body size.—In total, our dataset comprised 4,221 crick- 
ets, with > 1,100 crickets per sex for each morphological trait meas- 
ured, except for male tegmina color (Table 1). We first evaluated the 
relationship between three morphological traits that reflect overall 
body size in crickets: body length, femur length, and pronotum 
width. We found that body length measurements could vary de- 
pending on how crickets responded to being frozen in the lab or 
other factors such as number of eggs or last meal (see also Weissman 
and Gray 2019). Consequently, we chose to exclude body length 
measurements from our analyses but include them in our supple- 
mental datasets. Male and female individuals of both G. pennsylvan- 
icus and G. firmus had strong positive relationships between femur 
length and pronotum width (male G. pennsylvanicus: R? = 0.53, F, ,,, 
= 265, p < 2.2x10'° and male G. firmus: R’ = 0.76, F, ,, = 363.1, p 
< 2.2x10°'°, Suppl. material 1: fig. S1A; female G. pennsylvanicus: R? 
= 0.53, F, ,,, = 21, p < 2.2x10"° and female G. firmus: R’ = 0.74, F, ,, 
= 254.7, p < 2.2x10"°, Suppl. material 1: fig. S1B). Therefore, we 
used pronotum width as our estimate for overall body size to maxi- 
mize the number of individuals we could compare across datasets. 
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In female individuals, pronotum width and ovipositor length were 
also positively related in both species (G. pennsylvanicus: R* = 0.44, 
F, 4. = 165.7, p < 2.2x10"° and G. firmus: R* = 0.26, F, ,, = 30.48, 
p = 3.44x10’, Suppl. material 1: fig. S1C). In comparisons with 
G. thinos, we used femur length to estimate body size to maximize 


the number of individuals in those comparisons. 


Morphological differences between species.—There were significant 
differences among allopatric G. pennsylvanicus, G. firmus, G. thi- 
nos, and sympatric populations (e.g., G. firmus, G. pennsylvani- 
cus, and hybrids) in male body size (Kruskal-Wallis, 7? = 35.79, 
df = 3, p = 8.29x 10°), female body size (Kruskal-Wallis, 7? = 51.89, 
df = 3, p = 3.16x10"), female ovipositor length (Kruskal-Wallis, 
X = 1277.2, df = 3, p < 2.2x10"°), and relative ovipositor length 
(Kruskal-Wallis, 7? = 82.10, df = 3, p < 2.2x10"° ). When com- 
paring allopatric G. pennsylvanicus and G. firmus, male pronotum 
(p = 2.1x10°, Fig. 2A), female pronotum (p = 1.4x10", Fig. 2B), 
ovipositor length (p < 2.2x107°, Suppl. material 1: fig. S2A), and 
relative ovipositor length (p = 2.8x10°'°, Fig. 2C) were all signifi- 
cantly different. However, for each of these traits, there was still 
considerable overlap between allopatric species. Ovipositor length 
had the most striking differences between species (Suppl. material 
1: fig. S2A), even when controlling for body size (Fig. 2C). 

For males, tegmina color alone classified most individuals from 
allopatric populations as either G. pennsylvanicus or G. firmus (LDA, 
misclassification rate of 3%). One of the 24 G. pennsylvanicus males 
was misclassified as G. firmus, and zero of the 7 G. firmus males 
were misclassified as G. pennsylvanicus. When looking at male body 
size alone, the misclassification rate was much higher at 23%, with 
56 of the 268 G. pennsylvanicus males misclassified and 27 of the 
90 G. firmus males misclassified. There was not enough overlap in 
body size and tegmina color data to perform these analyses using 
both variables. For females, body size and relative ovipositor length 
classified most individuals from allopatric populations as either 
G. pennsylvanicus or G. firmus (LDA, misclassification rate 12%). 
Fifteen of the 189 G. pennsylvanicus were misclassified as G. firmus. 
and 17 of the 90 G. firmus were misclassified as G. pennsylvanicus. 


Table 1. Summary of sample sizes for morphological measurements by sex, population type, and region (See Fig. 1 for location information). 


Pronotum Width Femur Length Ovipositor Ovipositor Ovipositor Tegmina 
Females Males Females Males Length Pronotum Ratio Femur Ratio Color 
Totals 1203 1263 134 1213 4047 1174 1110 469 
CTR 4 5 4 5 12 4 4 - 
allopatric 4 5 4 5 4 4 4 - 
sympatric - - - 8 - - - 
NE 993 1010 871 849 3739 969 851 449 
allopatric 111 167 85 132 1480 108 82 23 
sympatric 882 843 786 717 2259 861 769 426 
NW 26 17 27 15 27 26 27 - 
allopatric 26 LZ 27 15 27 26 27, - 
SE 66 65 77 60 111 62 74 20 
allopatric 66 65 Oe 60 89 62 74 8 
sympatric - - - - 22 - - 18. 
SO 40 69 66 171 70 40 66 - 
allopatric 40 69 65 171 69 40 65 - 
sympatric - - 1 1 - 1 - 
SW 29 41 29 52 29 29 29 - 
allopatric 29 41 29 52 29 25 25 - 
WNC 45 56 60 61 59 44 59 - 
allopatric 45 56 60 61 59 44 ony - 
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Fig. 2. Allopatric populations of G. firmus and G. pennsylvanicus differ in overall body size and ovipositor length. A. Male pronotum 
width by species; B. Female pronotum width by species; C. Relative ovipositor length (ovipositor length/pronotum width). Boxplots 
indicate the mean values of each trait, quartiles, the range of the data (whiskers), and outliers. Individual data points are overlaid as 
scatterplots. Letters indicate the significant differences among groups (PWRST with corrected p-values < 0.05). 


Crickets from areas near the hybrid zone, which we refer to 
as sympatric, had considerable overlap with those from allopatric 
populations. Sympatric crickets were not different from G. firmus 
for male body size, but they were, on average, larger than G. pennsyl- 
vanicus (G. pennsylvanicus: p = 6.0x10°, G. firmus: p = 0.16, Fig. 2A) 
but were still different from both allopatric species for female 
body size (G. pennsylvanicus: p = 9.4x10”, G. firmus: p = 0.00032, 
Fig. 2B), female ovipositor length (G. pennsylvanicus: p < 2.0x10°1%, 
G. firmus: p < 2.0x10'°, Suppl. material 1: fig. S2A), and female 
relative ovipositor length (G. pennsylvanicus: p = 4.6x10°, G. fir- 
mus: p = 1.0x10°, Fig. 2C). This suggests that while these sympat- 
ric populations may have individuals that are more G. firmus-like 
or G. pennsylvanicus-like, they still have intermediate morphology 
compared to allopatric populations. 


Intraspecific variation in key morphological traits.—We then tested 
how these traits varied across the different geographic regions of 
each species. We found differences among regions of G. pennsyl- 
vanicus for male pronotum (Kruskal-Wallis, 7? = 56.11, df = 6, 
p = 2.7610"), female pronotum (Kruskal-Wallis, 7? = 63.44, df= 
6, p = 8.9x10°”), ovipositor length (Kruskal-Wallis, x? = 185.72, df 
= 6, p < 2.2x10"'°), and relative ovipositor length (Kruskal-Wallis, 
X = 33.6, df= 6, p = 8.03x10°%). Male and female G. pennsylvanicus 
were largest in the southern and midcentral US (SE, SO, SW, CTR, 
Fig. 3A, B), and they had the smallest body size in the northern 
west (WNC, NW). There were differences among regions in G. fir- 
mus pronotum width (Kruskal-Wallis, males, x? = 9.27, df = 2, 
p = 0.01; females, 7? = 9.15, df = 2, p = 0.01), in ovipositor length 
(Kruskal-Wallis, x? = 78.65, df = 2, p < 2.2x10'%), and relative 
ovipositor lengths (Kruskal-Wallis, 7? = 54.49, df = 2, p = 1.47x10- 
2), Male and female G. firmus were larger in the south than in the 
northeast, while G. firmus in the south were not significantly differ- 
ent from crickets in either the northeast or the southeast (Fig. 3A, 
B). In G. pennsylvanicus, ovipositor length varied by region. East- 
ern populations (NE, SE) had the shortest ovipositors, and those 
from the central US (CTR) had the longest ovipositors, although 
there was a very limited sample size for this region (Suppl. ma- 
terial 1: figs S2B, S3C). There was considerable variation in ovi- 
positor length among G. firmus populations; southern G. firmus 
females had significantly shorter relative ovipositors than G. fir- 
mus in the southeast, who in turn had shorter relative ovipositors 
than G. firmus in the northeast (Suppl. material 1: figs S2C, S3C). 
However, G. firmus in the southeast had very similar absolute ovi- 
positor lengths to northeastern G. firmus but had larger body sizes, 


whereas southern G. firmus simply had shorter ovipositors (Suppl. 
material 1: fig. S2C). 

Recent work by Weissman and Gray (2019) documented 
cryptic variation in southern USA G. firmus, so we took a closer 
look at these populations, separating crickets collected in Florida 
from those collected in Texas. We also included the recently de- 
scribed closely related species G. thinos, which is sympatric with 
Texas G. firmus (Weissman and Gray 2019). We found that male 
(Kruskal-Wallis, 7? = 29.26, df = 3, p = 1.98x10°) and female 
(Kruskal-Wallis, x* = 24.88, DF = 3, p = 1.63x10°) body size and 
ovipositor length (ovipositor length: Kruskal-Wallis, 7? = 101.39, 
df = 3, p< 2.2x10"%; relative ovipositor: Kruskal-Wallis, x? = 89.57, 
df = 3, p < 2.2x10"'°) differ among these groups (Fig. 4, Suppl. 
material 1: fig. S2). Compared to northeastern G. firmus, Florida 
G. firmus were much larger (Fig. 4A, B) but had only slightly larger 
ovipositor lengths (Suppl. material 1: fig. S2C), giving them short- 
er relative ovipositors (Fig. 4C). Texas G. firmus did not differ in 
overall body size from northeastern G. firmus but had even shorter 
relative ovipositor lengths (Fig. 4C, Suppl. material 1: fig. S2C). 
The magnitude of the morphological differences among Florida, 
Texas, and northeastern G. firmus is similar to that of the differenc- 
es between G. firmus and the recently described G. thinos. Gray et 
al. (2020) found that G. firmus in Texas and Florida are genetically 
distinct groups, with Texas G. firmus sister to G. pennsylvanicus and 
Florida G. firmus sister to both G. pennsylvanicus and Texas G. fir- 
mus. Altogether, the morphological differences and phylogenetic 
relationships support the findings by Weissman and Gray (2019) 
that Texas G. firmus may be an undescribed cryptic species. 


Morphology in sympatric populations.—For the subset of crickets that 
were from the hybrid zone or nearby (sympatric populations) and 
were also genotyped with molecular markers, we looked at the re- 
lationship between admixture and morphological traits. We found 
that each trait had a similar transition from G. pennsylvanicus to 
G. firmus, with highly admixed individuals having intermediate 
phenotypes (Fig. 5). We found that male pronotum (R?’ = 0.19, 
F, ,,) = 63.35, p = 4.3810"), male tegmina color (R? = 0.31, F, ,,, = 
60.82, p = 1.62x10"'”), female pronotum (R? = 0.28, F, ,,. = 107.3, 
p < 2.2x10"°), and relative ovipositor length (R’ = 0.47, F,,.. = 
243.1, p < 2.2x10"'%) all had strong correlation with the hybrid 
index. Because the SNPs used to calculate the hybrid index are con- 
centrated on the X chromosome (54 out of 110 (Maroja et al. 2015, 
Gainey et al. 2018)), females (XX) were more likely to be classi- 
fied with an intermediate hybrid index than males (XO). Overall, 
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Fig. 3. Cricket body size and relative ovipositor length varies by geographic region. A. Male pronotum width by species and region; 
B. Female pronotum width by species and region; C. Relative ovipositor length by species and region. Boxplots indicate the mean 
values of each trait, quartiles, the range of the data (whiskers), and outliers. Individual data points are overlaid as scatterplots. Letters 
indicate the significant differences among groups within each species (PWRST with corrected p-values < 0.05), and exact p-values are 
presented in Suppl. material 1: tables S1, S2. See Fig. 1 for location information. 


morphological traits were also correlated with mtDNA haplotypes: 
crickets that had G. pennsylvanicus mtDNA tended to be smaller 
(males: Kruskal-Wallis, x? = 43.14, df= 1, p = 5.1110"; females: 
Kruskal-Wallis, x* = 44.86, df= 1, p = 2.11x10-"), darker (Kruskal- 
Wallis, x? = 33.75, df = 1, p = 6.27x10°) crickets with shorter ovi- 
positors (Kruskal-Wallis, x* = 37.67, df= 1, p = 8.40x 10°’) (Fig. 6). 
We found that crickets with G. firmus ancestry at nuclear markers 
(hybrid index = 1) often had G. pennsylvanicus mtDNA haplotypes 
(Fig. 7), indicating asymmetric introgression of the mtDNA. 


Environmental predictors of morphology.—In allopatric populations 
throughout broad ranges, we found that latitude, elevation, average 


soil percent clay, and minimum and maximum temperatures cre- 
ated the best model for ovipositor length. Latitude, longitude, soil 
percent sand, and minimum temperature created the best model 
for pronotum width (Table 2). Average soil percent clay and higher 
minimum and maximum temperatures were positively associated 
with longer ovipositor lengths, and higher minimum temperatures 
were positively associated with larger body size, which are char- 
acteristics of G. firmus (Suppl. material 1: fig. $3). In the subset 
of Connecticut sympatric and allopatric populations, minimum 
and maximum temperatures, as well as soil percent organic mat- 
ter, created the best model, with positive associations for all three 
variables and ovipositor length (Table 2, Suppl. material 1: fig. $3). 
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Fig. 4. Morphological variation in G. firmus consistent with proposed cryptic species. A. Male femur length; B. Female femur length; 
C. Relative ovipositor length (ovipositor length/femur length). There is considerable morphological variation among northeastern, 
Florida, and Texas G. firmus, which is similar to the magnitude of morphological divergence observed in the closely related species 
G. thinos. This combined with genetic divergence suggests there may be cryptic species in what is currently considered G. firmus. Box- 
plots indicate the mean values of each trait, quartiles, the range of the data (whiskers), and outliers. Individual data points are overlaid 
as scatterplots. Letters indicate the significant differences among groups (PWRST with corrected p-values < 0.05), and exact p-values are 
presented in Suppl. material 1: table S3. 
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Fig. 5. Crickets with more hybrid background have intermediate morphological traits. The relationship between the hybrid index (an 
estimate of ancestry proportions, G. pennsylvanicus = 0 and G. firmus = 1) and A. Male pronotum width; B. Female pronotum width; 
C. Relative ovipositor length, and D. Male tegmina color. 
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Fig. 6. Morphological traits tended to correspond to mtDNA haplotypes. A. Male pronotum width; B. Female pronotum width; C. Rela- 
tive ovipositor length; and D. Male tegmina color. Boxplots indicate the mean values of each trait, quartiles, the range of the data 
(whiskers), and outliers. Individual data points are overlaid as scatterplots. 


Table 2. Results of linear regression and AIC to test the relationship between environmental variables and morphological traits in 
female crickets of both species. | indicates variables where values are based on the year the samples were collected. 


A. Ovipositor length 


Df Sum of Sq RSS AIC Coefficient St. Error t-value p-value 
(Intercept) - 887.86 321.96 16.213 0.132 122.417 < 2.00E-16 
Latitude 1 9.283 897.14 322.33 0.545 0.358 1523 0.129 
Precipitation!’ 1 1.054 886.81 323.69 - - - - 
Longitude 1 0.389 887.47 323.86 - - - - 
Human Footprint 1 0.132 887.73 323.93 - - - - 
Avg Soil % Sand 1 0.015 887.85 323.96 - - - - 
Avg Soil % Organic Matter 1 0.004 887.86 323,96 - - - - 
Elevation 1 26.035 913.90 326.55 -0.638 0.250 -2.551 0.011 
Avg Soil % Clay 1 26.562 914.42 326.68 0.365 0.142 2.577 0.011 
Minimum Temperature’ 1 29.311 SUASLT, 327.36 -1.244 0.459 -2.707 0.007 
Maximum Temperature’ 1 124.629 1012.49 349.91 2.281 0.409 5.582 6.89E-08 

B. Pronotum width 

(Intercept) - 56.539 25309 5.83503 0.036 162.636 < 2.00E-16 
Maximum Temperature’ 1 0.381 36.158 -253.69 - - - - 
Precipitation! 1 0.176 36.363 20 2;¢3 - - - - 
Human Footprint 1 0.147 36.392 -252.59 - - - - 
Elevation 1 0.103 36.436 -252.38 - - - - 
Avg Soil % Clay 1 0.088 36.451 -252.31 - - - - 
Avg Soil % Organic Matter 1 0.013 36.526 -251.96 - - - - 
Minimum Temperature’ 1 2.964 39.503 -242.56 -0.233 0.064 -3.670 3.27E-04 
Avg Soil % Sand 1 3953 40.492 -238.33 -0.180 0.043 -4.238 3.73E-05 
Longitude 1 4.618 41.157 -235.55 -0.187 0.041 -4.580 9.07E-06 
Latitude 1 12.890 49.429 -204.23 -0.536 0.070 -7.652 1.53E-12 
Discussion 2009a), can be largely differentiated by a combination of body 


Cryptic diversity in a wide-ranging species. —The hybrid zone between 
the field crickets G. firmus and G. pennsylvanicus has been a model 
for understanding speciation (Harrison and Rand 1989, Harrison 
and Larson 2014). The field cricket hybrid zone stretches from the 
northeastern United States as far south as Virginia and likely far- 
ther into the southeast. Divergence in morphology, nuclear and 
mitochondrial DNA, and reproductive barriers have been careful- 
ly studied in several major regions of the hybrid zone (Harrison 
1985, Rand and Harrison 1989, Ross and Harrison 2002, Maroja 
et al. 2009a, 2009b, Larson et al. 2012, 2014). Yet, even in this 
well-studied system, there is geographic diversity across the ranges 
of these species that complicates their relationships. 

Our results confirm that allopatric populations of these two 
species, defined by genetic markers (Harrison and Arnold 1982, 
Willett et al. 1997, Broughton and Harrison 2003, Maroja et al. 


size, male tegmina color, and female ovipositor length (Fig. 2). 
At the same time, there is regional variation in these traits with- 
in each species (Fig. 3). These differences may be due to local 
adaptation of life history traits such as egg diapause and devel- 
opment time (discussed below) or phenotypic plasticity. How- 
ever, in some cases, they may also indicate cryptic diversity in 
field crickets. 

In their revision of North American field crickets, Weissman 
and Gray (2019) proposed that there was cryptic diversity in the 
southern populations of G. firmus, particularly in Texas. Impor- 
tantly, our phenotypic comparisons confirmed that Texas and 
Florida G. firmus are morphologically distinct from northeastern 
G. firmus (Fig. 4). In a recent nuclear phylogeny, Texas and Florida 
G. firmus-like crickets also formed distinct clusters within the larg- 
er G. pennsylvanicus group (Weissman and Gray 2019, Gray et al. 
2020). Unfortunately, we do not have a phylogeny that includes 
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Fig. 7. Mitochondrial DNA introgression is largely asymmetric. 
Crickets with G. firmus ancestry at nuclear markers (hybrid in- 
dex = 1) often had G. pennsylvanicus mtDNA haplotypes. 


genes from both Texas and Florida G. firmus and northeastern 
G. firmus, so the relationships among these groups are still un- 
clear. However, the combination of distinct morphology and phy- 
logenetic relationships suggests that at least one cryptic species 
of Gryllus exists, a situation that will not be resolved without fur- 
ther genotyping and/or evaluations of reproductive compatibility 
among these populations. 


Intermediate phenotypes in hybrid zone crickets. —The morphologi- 
cal traits that best distinguish species in allopatry can also be 
used to distinguish these species in or near the hybrid zone. In 
this study, we took a conservative approach to defining allopat- 
ric and sympatric populations. Allopatric populations were 
those well outside of where the two species co-occur and are 
typically populations that have been genotyped with species- 
diagnostic markers. We found that in sympatry, crickets that 
were mostly G. firmus or mostly G. pennsylvanicus at nuclear 
markers (Larson et al. 2013a, 2014) had morphological traits 
that are also G. firmus-like or G. pennsylvanicus-like. Both male 
and female body size, male tegmina color, and relative oviposi- 
tor length had clinal variations from G. pennsylvanicus-like to 
G. firmus-like, with highly admixed individuals having inter- 
mediate phenotypes (Fig. 5). Male tegmina color stood out as 
having the fewest individuals with intermediate hybrid index 
values (Fig. 5D), but this is most likely because the SNPs used 
to calculate the hybrid index were predominately X-linked, so 
male XO crickets were rarely heterozygous at those SNPs and 
had overall lower hybrid indices (Larson et al. 2014, Maroja et 
al. 2015, Gainey et al. 2018). 

The relationship between morphology and mitochondrial 
haplotype was less clear for populations near or in the hybrid 
zone. Crickets that were mostly G. firmus at the nuclear mark- 
ers often had G. pennsylvanicus mtDNA (Fig. 7). This pattern 
fits with what we expect based on the one-way prezygotic in- 
compatibility between G. firmus females and G. pennsylvanicus 
males (Harrison 1983, Maroja et al. 2009b, Larson et al. 2012). 
All F1 hybrids are produced from crosses with G. pennsylvanicus 
mothers; thus, G. pennsylvanicus mtDNA will be more likely to 
introgress into G. firmus. Even rare instances of hybridization 
might lead to mtDNA introgression, such as the mtDNA capture 
observed in many mammal species (Melo-Ferreira et al. 2005, 
Good et al. 2008). 
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Adaptations to soil type.—Ovipositor length is one of the most 
striking morphological differences between G. firmus and 
G. pennsylvanicus. Female crickets use their ovipositors to lay their 
eggs in the soil, and ovipositor length has been hypothesized to re- 
late to the soil type and/or the depth of egg laying (Masaki 1979). 
The depth of egg laying may be a particularly critical life-history 
trait in G. pennsylvanicus and G. firmus because these species over- 
winter as eggs, as opposed to most field crickets that overwinter 
as early instar nymphs (Alexander 1968, Harrison and Bogdano- 
wicz 1995). For eggs to be viable, they must withstand low winter 
temperatures and freeze/thaw cycles (Ross and Harrison 2006). 
Throughout its range, G. firmus is most often found on sandy 
coastal soils (Harrison and Arnold 1982, Weissman and Gray 
2019) and tends to have a longer ovipositor than G. pennsylvanicus 
(Fig. 3, Suppl. material 1: fig. S2). This may be an adaptation to 
laying eggs deeper in sandy substrates in response to intermittent 
rainfall and the risk of eggs drying out (Walker 1980). In some 
parts of the hybrid zone, such as Connecticut, the association with 
different soil types is striking. The two species have been found 
on micro habitat patches of loam (G. pennsylvanicus) and sandy 
(G. firmus) soils in Connecticut (Harrison 1986, Harrison and 
Rand 1989, Rand and Harrison 1989), and interactions between 
the two species occur across these habitat patch boundaries on a 
scale of only hundreds of meters (Ross and Harrison 2002, Larson 
et al. 2014). 

Despite what appears to be strong habitat associations, the 
relationship between soil type and ovipositor length is compli- 
cated. Ovipositor length does not necessarily determine egg-laying 
depth; instead, females may wield long ovipositors at different an- 
gles (Réale and Roff 2002). It is also not clear exactly how the as- 
sociation between ovipositor length and soil type is maintained. 
Females of both species prefer to lay eggs in loamy soil, and there 
is no difference in overwintering egg viability in different soil 
types (Ross and Harrison 2006). Finally, these associations are 
clearly established only in a small part of the species’ ranges, i.e., 
Connecticut (Rand and Harrison 1989, Ross and Harrison 2002, 
Larson et al. 2013b). Even where soil associations appear to be the 
strongest, the transition from sandy to loamy soils is more gradual 
and less distinct than we might expect based on the patchiness 
of G. firmus and G. pennsylvanicus populations (Ross and Harri- 
son 2002, Larson et al. 2014). Here we find that both across the 
broad ranges of these species and at an intermediate scale in the 
Connecticut hybrid zone, there is no strong association between 
ovipositor length and sandy soils. In fact, we tend to see crickets 
with longer ovipositors on clay soils (Table 2). This might be due 
to the different methods used to quantify soil type (soil survey 
data versus on-site soil sampling), but altogether, this suggests that 
habitat associations in these species are variable and should be 
investigated further. 


Body size, climate, and life cycle.—In insects, seasonality and the 
length of the growing season are critical to the rate of development 
and adult body size (Masaki 1961, Tauber and Tauber 1981). This 
is particularly true for hemimetabolous insects, which often go 
through many nymphal stages and have long development times 
before reaching their full size and sexual maturity (Kivela et al. 
2011). Insects at higher latitudes have shorter growing seasons and, 
as a result, may develop more quickly or reach an overall smaller 
body size (Masaki 1967, Parsons and Joern 2014). This pattern 
of smaller body sizes at higher latitudes is sometimes referred to 
as the converse of Bergman’s rule, which states that individuals 
have larger body sizes in colder climes (Masaki 1967, Mousseau 
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1997). We see this pattern most clearly in G. pennsylvanicus, where 
we found that populations with the smallest body sizes tended to 
be farther north (WNC and NW, Fig. 2). Indeed, we found that 
crickets at higher latitudes had, on average, smaller body sizes and 
that there was a significant relationship between body size and 
latitude (Table 2). 

We may not expect a direct relationship between body size 
and latitude if the length of the growing season allows for mul- 
tiple generations per year. Insects can shift from continuous de- 
velopment in the south to univoltine (one generation per year) in 
the north (Masaki 1961, 1967). As a result, there may be regions 
where body size is smaller than expected based on latitude to ac- 
commodate multiple generations per year. We did not find this 
pattern in our results, but we may not have had the resolution of 
latitudinal samples to see a sawtooth pattern in body size. Howev- 
er, there is some evidence that development time in G. firmus var- 
ies with latitude. In Virginia, G. firmus emerge earlier in the season 
than G. pennsylvanicus, leading to temporal isolation in that part 
of the hybrid zone, but in Connecticut, the two emerge simultane- 
ously (Harrison 1985). In Florida, G. firmus is reported to have 
multiple generations per year (Walker, personal observation; re- 
ported in Weissman and Gray 2019), where throughout its range, 
it otherwise appears to have a single generation per year (Walker 
1980). Notably, despite having many generations per year, Florida 
G. firmus are considerably larger than northern populations. It is 
unclear whether there is a continuous shift in life cycle across the 
range of G. firmus or if Florida G. firmus have a distinct life history 
from other G. firmus. 


Conclusions 


In studies of speciation and to understand the effects of local 
selection, it is critical to quantify morphological and genetic varia- 
tions across the geographic range of widespread species. The field 
cricket hybrid zone is an example of how important the larger geo- 
graphic context can be. In some regions of the field cricket hybrid 
zone, G. pennsylvanicus and G. firmus have a patchy distribution, 
and G. firmus crickets are found on sandy soils (Rand and Har- 
rison 1989, Ross and Harrison 2002). However, the strong soil as- 
sociation breaks down in other regions of the hybrid zone (Larson 
et al. 2013b) and across their geographic range, suggesting that the 
soil association may be a result of local adaptation or colonization 
history (Hauffe and Searle 1993, Gompert et al. 2010). Our results 
provide a foundation for future geographically expansive studies 
that compare genetic divergence and the role of specific traits in 
reproductive barriers to better understand local adaptation and 
speciation in this system. More broadly, this is an example of how 
critical it is to move studies of speciation beyond the comparison 
of a few focal populations. Geographically expansive studies of 
phenotypic and genetic divergence will also be important for un- 
derstanding how species distributions and hybrid zones shift over 
time and in a changing climate (Britch and Cain 2001, Taylor et 
al. 2015). 


Author Contributions 


ARB, CJ, DBW, DAG, CLR, LSM and ELL collected the data; 
ARB, CJ, and ELL combined the datasets; AG obtained the envi- 
ronmental data and advised on analyses. ARB conducted all analy- 
ses. ARB and ELL wrote the manuscript with contributions from 
all authors. 


A.R. BYERLY, C. JENCK, A.R.B. GOETZ, D.B. WEISSMAN, D.A. GRAY, C.L. ROSS, L.S. MAROJA AND E.L. LARSON 


Acknowledgements 


Rick Harrison collected many of the samples in this dataset. 
He was a wonderful mentor and friend to many of the authors 
and an inspiration to all of us. We thank the University of Den- 
ver Ecology and Evolution group, especially Shannon Murphy 
and Jonathan Velotta, for their feedback on this manuscript. This 
work was supported by grants to ARB from the Explorer's Club, 
The Orthopterists’ Society, and the Society of Systematic Biologists 
and National Science Foundation grants to ELL (DEB 2012041) 
and LSM (DEB 2012060). 


References 


Alexander RD (1968) Life cycle origins, speciation, and related phenom- 
ena in crickets. The Quarterly Review of Biology 43: 1-41. https://doi. 
org/10.1086/405628 

Andrés JA, Larson EL, Bogdanowicz SM, Harrison RG (2013) Patterns of 
transcriptome divergence in the male accessory gland of two closely 
related species of field crickets. Genetics 193: 501-513. https://doi. 
org/10.1534/genetics.112.142299 

Barton NH, Hewitt GM (1985) Analysis of hybrid zones. Annual Review 
of Ecology and Systematics 16: 113-148. https://doi.org/10.1146/an- 
nurev.es. 16.110185.000553 

Benjamini Y, Hochberg Y (1995) Controlling the false discovery rate: a 
practical and powerful approach to multiple testing. Journal of 
the Royal Statistical Society. Series B, Statistical Methodology 57: 
289-300. https://doi.org/10.1111/j.2517-6161.1995.tb02031.x 

Britch SC, Cain ML (2001) Spatio-temporal dynamics of the Allonemo- 
bius fasciatus— A. socius mosaic hybrid zone: a 14-year perspective. 
Molecular Ecology 10: 627-638. https://doi.org/10.1046/j.1365- 
294x.2001.01215.x 

Broughton RE, Harrison RG (2003) Nuclear gene genealogies reveal his- 
torical, demographic and selective factors associated with speciation 
in field crickets. Genetics 163: 1389-1401. https://doi.org/10.1093/ 
genetics/163.4.1389 

Gainey DP, Kim JY, Maroja LS (2018) Mapping reduced introgression 
loci to the X chromosome of the hybridizing field crickets, Gryllus 
firmus and G. pennsylvanicus. PLoS ONE 13: e0208498. https://doi. 
org/10.1371/journal.pone.0208498 

Gavrilets S, Arnqvist G, Friberg U (2001) The evolution of female mate 
choice by sexual conflict. Proceedings of the Royal Society of London. 
Series B: Biological Sciences 268: 531-539. https://doi.org/10.1098/ 
rspb.2000.1382 

Gompert Z, Lucas LK, Fordyce JA, Forister ML, Nice CC (2010) Sec- 
ondary contact between Lycaeides idas and L. melissa in the Rocky 
Mountains: extensive admixture and a patchy hybrid zone. Mo- 
lecular Ecology 19: 3171-3192. https://doi.org/10.1111/j.1365- 
294X.2010.04727.x 

Gonzalez C, Ornelas JE Gutiérrez-Rodriguez C (2011) Selection and 
geographic isolation influence hummingbird speciation: genetic, 
acoustic and morphological divergence in the wedge-tailed sabrew- 
ing (Campylopterus curvipennis). BMC Evolutionary Biology 11: 38. 
https://doi.org/10.1186/1471-2148-11-38 

Good JM, Hird S, Reid N, Demboski JR, Steppan SJ, Martin-Nims TR, Sul- 
livan J (2008) Ancient hybridization and mitochondrial capture be- 
tween two species of chipmunks. Molecular Ecology 17: 1313-1327. 
https://doi.org/10.1111/j.1365-294X.2007.03640.x 

Gray DA, Weissman DB, Cole JA, Lemmon EM, Lemmon AR (2020) Mul- 
tilocus phylogeny of Gryllus field crickets (Orthoptera: Gryllidae: 
Gryllinae) utilizing anchored hybrid enrichment. Zootaxa 4750: 
328-348. https://doi.org/10.11646/zootaxa.4750.3.2 

Harrison RG (1983) Barriers to gene exchange between closely relat- 
ed cricket species. I. Laboratory hybridization studies. Evolution; 
International Journal of Organic Evolution 37: 245-251. https://doi. 
org/10.2307/2408333 


JOURNAL OF ORTHOPTERA RESEARCH 2023, 32(2) 


A.R. BYERLY, C. JENCK, A.R.B. GOETZ, D.B. WEISSMAN, D.A. GRAY, C.L. ROSS, L.S. MAROJA AND E.L. LARSON 


Harrison RG (1985) Barriers to gene exchange between closely related 
cricket species. II. Life cycle variation and temporal isolation. Evo- 
lution; International Journal of Organic Evolution 39: 244-259. 
https://doi.org/10.2307/2408360 

Harrison RG (1986) Pattern and process in a narrow hybrid zone. Heredity 
56: 337-349. https://doi.org/10.1038/hdy.1986.55 

Harrison RG (1990) Hybrid zones: windows on evolutionary process. 
Oxford Surveys in Evolutionary Biology 7: 69-128. 

Harrison RG, Arnold J (1982) A narrow hybrid zone between closely related 
cricket species. Evolution 36: 535. https://doi.org/10.2307/2408099 

Harrison RG, Rand DM (1989) Mosaic hybrid zones and the nature of spe- 
cies boundaries. Speciation and its Consequences: 111-133. 

Harrison RG, Bogdanowicz SM (1995) Mitochondrial DNA phylogeny of 
North American field crickets: perspectives on the evolution of life cy- 
cles, songs, and habitat associations. Journal of Evolutionary Biology 
8: 209-232. https://doi.org/10.1046/j.1420-9101.1995.8020209.x 

Harrison RG, Bogdanowicz SM (1997) Patterns of variation and link- 
age disequilibrium in a field cricket hybrid zone. Evolution 51: 493. 
https://doi.org/10.2307/2411122 

Harrison RG, Larson EL (2014) Hybridization, introgression, and the na- 
ture of species boundaries. The Journal of Heredity 105: 795-809. 
https://doi.org/10.1093/jhered/esu033 

Harrison RG, Larson EL (2016) Heterogeneous genome divergence, dif- 
ferential introgression, and the origin and structure of hybrid zones. 
Molecular Ecology 25: 2454-2466. https://doi.org/10.1111/mec.13582 

Harrison RG, Rand DM, Wheeler WC (1987) Mitochondrial DNA varia- 
tion in field crickets across a narrow hybrid zone. Molecular Biology 
and Evolution 4: 144-144. 

Hauffe HC, Searle JB (1993) Extreme karyotypic variation in a Mus muscu- 
lus domesticus hybrid zone: the tobacco mouse story revisited. Evo- 
lution; International Journal of Organic Evolution 47: 1374-1395. 
https://doi.org/10.2307/2410154 

Hollander J, Montano-Rend6n M, Bianco G, Yang X, Westram AM, Duvaux 
L, Reid DG, Butlin RK (2018) Are assortative mating and genital di- 
vergence driven by reinforcement? Evolution Letters 2: 557-566. 
https://doi.org/10.1002/ev13.85 

Jenck CS, Lehto WR, Ketterman BT, Sloan LE, Sexton AN, Tinghitella RM 
(2020) Phenotypic divergence among threespine stickleback that 
differ in nuptial coloration. Ecology and Evolution 10: 2900-2916. 
https://doi.org/10.1002/ece3.6105 

Jiménez RA, Ornelas JF (2015) Historical and current introgression in a 
Mesoamerican hummingbird species complex: a biogeographic per- 
spective. Peer) 4: e1556. https://doi.org/10.7717/peerj.1556 

Karl T, Koss WJ (1984) Regional and national monthly, seasonal, and an- 
nual temperature weighted by area, 1895-1983. https://repository. 
library.noaa.gov/view/noaa/ 10238 [April 12, 2021] 

Kivela SM, Valimaki P, Carrasco D, Maenpéa MI, Oksanen J (2011) Lati- 
tudinal insect body size clines revisited: a critical evaluation of the 
saw-tooth model. The Journal of Animal Ecology 80: 1184-1195. 
https://doi.org/10.1111/j.1365-2656.2011.01864.x 

Larson EL, Hume GL, Andrés JA, Harrison RG (2012) Post-mating prezy- 
gotic barriers to gene exchange between hybridizing field crickets. 
Journal of Evolutionary Biology 25: 174-186. https://doi.org/10.1111/ 
j.1420-9101.2011.02415.x 

Larson EL, Andrés JA, Bogdanowicz SM, Harrison RG (2013a) Differential 
introgression in a mosaic hybrid zone reveals candidate barrier genes. 
Evolution; International Journal of Organic Evolution 67: 3653- 
3661. https://doi.org/10.1111/evo.12205 

Larson EL, Guilherme Becker C, Bondra ER, Harrison RG (2013b) Structure 
of a mosaic hybrid zone between the field crickets Gryllus firmus and 
G. pennsylvanicus. Ecology and Evolution 3: 985-1002. https://doi. 
org/10.1002/ece3.514 

Larson EL, White TA, Ross CL, Harrison RG (2014) Gene flow and the 
maintenance of species boundaries. Molecular Ecology 23: 1668- 
1678. https://doi.org/10.1111/mec.12601 

Maroja LS, Andrés JA, Harrison RG (2009a) Genealogical discordance 
and patterns of introgression and selection across a cricket hybrid 


199 


zone. Evolution; International Journal of Organic Evolution 63: 
2999-3015. https://doi.org/10.1111/j.1558-5646.2009.00767.x 

Maroja LS, Andrés JA, Walters JR, Harrison RG (2009b) Multiple barriers 
to gene exchange in a field cricket hybrid zone. Biological Journal 
of the Linnean Society 97: 390-402. https://doi.org/10.1111/j.1095- 
8312.2009.01201.x 

Maroja LS, Larson EL, Bogdanowicz SM, Harrison RG (2015) Genes with 
restricted introgression in a field cricket (Gryllus firmus/Gryllus penn- 
sylvanicus) hybrid zone are concentrated on the X chromosome and 
a single autosome. G3 Genes|Genomes|Genetics 5: 2219-2227. 
https://doi.org/10.1534/g3.115.021246 

Masaki S (1961) Geographic variation of diapause in insects. Bulletin of 
the Faculty of Agriculture, Hirosaki University 7: 66-98. 

Masaki S (1967) Geographic variation and climatic adaptation in a field 
cricket (Orthoptera: Gryllidae). Evolution; International Journal of 
Organic Evolution 21: 725-741. https://doi.org/10.2307/2406770 

Masaki S (1979) Climatic adaptation and species status in the lawn ground 
cricket. Oecologia 43: 207-219. https://doi.org/10.1007/BF00344771 

Melo-Ferreira J, Boursot P, Suchentrunk FE, Ferrand N, Alves PC (2005) 
Invasion from the cold past: extensive introgression of mountain 
hare (Lepus timidus) mitochondrial DNA into three other hare spe- 
cies in northern Iberia. Molecular Ecology 14: 2459-2464. https:// 
doi.org/10.1111/j.1365-294X.2005.02599.x 

Moran PA, Hunt J, Mitchell C, Ritchie MG, Bailey NW (2020) Sexual se- 
lection and population divergence III: Interspecific and intraspecific 
variation in mating signals. Journal of Evolutionary Biology 33: 990- 
1005. https://doi.org/10.1111/jeb.13631 

Mousseau TA (1997) Ecotherms follow the converse to Bergmann’s rule. 
Evolution; International Journal of Organic Evolution 51: 630-632. 
https://doi.org/10.2307/2411138 

Oneal E, Knowles LL (2013) Ecological selection as the cause and sexual 
differentiation as the consequence of species divergence? Proceedings 
of the Royal Society B, Biological Sciences 280: 20122236. https://doi. 
org/10.1098/rspb.2012.2236 

Parsons SMA, Joern A (2014) Life history traits associated with body size 
covary along a latitudinal gradient in a generalist grasshopper. Oeco- 
logia 174: 379-391. https://doi.org/10.1007/s00442-013-2785-6 

Poggio L, de Sousa LM, Batjes NH, Heuvelink GBM, Kempen B, Ribeiro E, 
Rossiter D (2021) SoilGrids 2.0: producing soil information for the 
globe with quantified spatial uncertainty. Soil 7: 217-240. https://doi. 
org/10.5194/soil-7-217-2021 

Polly PD, Wojcik JM (2019) Geometric morphometric tests for phenotypic 
divergence between chromosomal races. In: Shrews, Chromosomes 
and Speciation. Cambridge University Press, 336-364. https://doi. 
org/10.1017/9780511895531.011 

Rand DM, Harrison RG (1989) Ecological genetics of a mosaic hybrid 
zone: mitochondrial, nuclear and reproductive differentiation of 
crickets by soil type. Evolution; International Journal of Organic Evo- 
lution 43: 432-449. https://doi.org/10.2307/2409218 

R Core Team (2020) R: A language and environment for statistical com- 
puting. R Foundation for Statistical Computing, Vienna, Austria. 
https://www.R-project.org/ 

Réale D, Roff DA (2002) Quantitative genetics of oviposition behaviour 
and interactions among oviposition traits in the sand cricket. Animal 
Behaviour 64: 397-406. https://doi.org/10.1006/anbe.2002.3084 

Ross CL, Harrison RG (2002) A fine-scale spatial analysis of the mosaic 
hybrid zone between Gryllus firmus and Gryllus pennsylvanicus. Evo- 
lution; International Journal of Organic Evolution 56: 2296-2312. 
https://doi.org/10.1111/j.0014-3820.2002.tb00153.x 

Ross CL, Harrison RG (2006) Viability selection on overwintering eggs 
in a field cricket mosaic hybrid zone. Oikos 115: 53-68. https://doi. 
org/10.1111/j.2006.0030-1299.15054.x 

Runemark A, Hansson B, Pafilis P, Valakos ED, Svensson EI (2010) Is- 
land biology and morphological divergence of the Skyros wall liz- 
ard Podarcis gaigeae: a combined role for local selection and ge- 
netic drift on color morph frequency divergence? BMC Evolutionary 
Biology 10: 269. 


JOURNAL OF ORTHOPTERA RESEARCH 2023, 32(2) 


200 


Shastry V, Adams PE, Lindtke D, Mandeville EG, Parchman TL, Gompert 
Z, Buerkle CA (2021) Model-based genotype and ancestry estimation 
for potential hybrids with mixed-ploidy. Molecular Ecology Resources 
21: 1434-1451. https://doi.org/10.1111/1755-0998.13330 

Shaw KL, Mullen SP (2011) Genes versus phenotypes in the study of spe- 
ciation. Genetica 139: 649-661. https://doi.org/10.1007/s10709-011- 
9562-4 

Skoglund S, Siwertsson A, Amundsen P-A, Knudsen R (2015) Morphologi- 
cal divergence between three Arctic charr morphs - the significance of 
the deep-water environment. Ecology and Evolution 5: 3114-3129. 
https://doi.org/10.1002/ece3.1573 

Sottas C, Reif J, Kuczynski L, Reifova R (2018) Interspecific competition 
promotes habitat and morphological divergence in a secondary con- 
tact zone between two hybridizing songbirds. Journal of Evolutionary 
Biology 31: 914-923. https://doi.org/10.1111/jeb.13275 

Tauber CA, Tauber MJ (1981) Insect seasonal cycles. Annual Review of 
Ecology and Systematics 12: 281-308. https://doi.org/10.1146/an- 
nurev.es. 12.110181.001433 

Taylor SA, Larson EL, Harrison RG (2015) Hybrid zones: windows on cli- 
mate change. Trends in Ecology & Evolution 30: 398-406. https://doi. 
org/10.1016/j.tree.2015.04.010 

Uyeda JC, Arnold SJ, Hohenlohe PA, Mead LS (2009) Drift promotes 
speciation by sexual selection. Evolution; International Journal of 
Organic Evolution 63: 583-594. https://doi.org/10.1111/j.1558- 
5646.2008.00589.x 

Venter O, Sanderson EW, Magrach A, Allan JR, Beher J, Jones KR, Possin- 
gham HP, Laurance WE, Wood P, Fekete BM, Levy MA, Watson JEM 
(2016) Sixteen years of change in the global terrestrial human foot- 
print and implications for biodiversity conservation. Nature Commu- 
nications 7: 12558. https://doi.org/10.1038/ncomms12558 

Venter O, Sanderson EW, Magrach A, Allan JR, Beher J, Jones KR, Poss- 
ingham HP, Laurance WE, Wood P, Fekete BM, Levy MA, Watson JE 
(2018) Last of the Wild Project, Version 3 (LWP-3): 2009 Human 
Footprint, 2018 Release. Palisades, New York: NASA Socioeconomic 
Data and Applications Center (SEDAC). https://doi.org/10.7927/ 
H46TOJQ4 [Accessed October 18, 2021] 

Walker TJ (1980) Mixed oviposition in individual females of Gryllus firmus: 
Graded proportions of fast-developing and diapause eggs. Oecologia 
47: 291-298. https://doi.org/10.1007/BF00398519 

Wang G-H, Li H, Zhao H-W, Zhang W-K (2017) Detecting climatically 
driven phylogenetic and morphological divergence among spruce 
(Picea) species worldwide. Biogeosciences 14: 2307-2319. https:// 
doi.org/10.5194/bg-14-2307-2017 

Weissman DB, Gray DA (2019) Crickets of the genus Gryllus in the United 
States (Orthoptera: Gryllidae: Gryllinae). Zootaxa 4705: 1- 277. htt- 
ps://doi.org/10.11646/zootaxa.4705.1.1 

Willett CS, Ford MJ, Harrison RG (1997) Inferences about the origin of a 
field cricket hybrid zone from a mitochondrial DNA phylogeny. He- 
redity 79: 484-494. https://doi.org/10.1038/hdy.1997.188 

Wolf JBW, Harrod C, Brunner S, Salazar S, Trillmich F, Tautz D (2008) 
Tracing early stages of species differentiation: ecological, morpho- 
logical and genetic divergence of Galapagos sea lion populations. 
BMC Evolutionary Biology 8: 150. https://doi.org/10.1186/1471- 
2148-8-150 


A.R. BYERLY, C. JENCK, A.R.B. GOETZ, D.B. WEISSMAN, D.A. GRAY, C.L. ROSS, L.S. MAROJA AND E.L. LARSON 


Supplementary material 1 


Author: Amy R. Byerly, Clara Jenck, Alexander R. B. Goetz, Da- 

vid B. Weissman, David A. Gray, Charles L. Ross, Luana S. Maroja, 

Erica L. Larson 

Data type: docx 

Explanation note: table $1. P-values for PWRST posthoc contrasts of 
allopatric G. pennsylvanicus populations by region (See Fig. 3). P- 
values marked with *** lost significance after correction. table 
S2. P-values for PWRST posthoc contrasts of allopatric G. firmus 
populations by region (See Fig. 3). P-values marked with *** lost 
significance after correction. table $3. P-values for PWRST post- 
hoc contrasts of G. thinos populations and G. firmus populations 
in the Northeast, Florida, and Texas. See Fig. 4. figure S1. Rela- 
tionship among phenotypic characteristics for allopatric popu- 
lations. (a) Male Gryllus firmus (yellow) and G. pennsylvanicus 
(teal) linear models show positive relationships in femur length 
and pronotum width (G. firmus: R2 = 0.76, F1,117 = 363.1, p-val- 
ue < 2.2e-16; G. pennsylvanicus: R2 = 0.53, F1,233 = 265.0, p-value 
< 2.2e-16). Female linear models show positive relationships in 
both (b) femur length and pronotum width (G. firmus: R2 = 0.74, 
F1,89 = 254.7,p-value < 2.2e-16; G. pennsylvanicus: R2 = 0.53, 
F1,192 = 217.0, p-value < 2.2e-16) and (c) ovipositor length and 
pronotum width (G. firmus: R2 = 0.26, F1,87 = 30.48, p-value = 
3.44e-07; G. pennsylvanicus: R2 = 0.44, F1,214 = 165.7, p-value < 
2.2e-16). figure S2. Ovipositor length differences between species 
and among populations of each species. A. Ovipositor length dif- 
ferences between G. firmus, G. pennsylvanicus and sympatric popu- 
lations (G. firmus vs G. pennsylvanicus p < 2.0E-16; G. firmus vs 
sympatric p <2.0E-16; G. pennsylvanicus vs sympatric p <2.0E-16). 
B. Ovipositor length differences among populations of G. penn- 
sylvanicus. Posthoc p-values are presented in Suppl. material 1: 
table S1. C. Ovipositor length differences among populations of 
G. firmus and G. thinos. Posthoc p-values are presented in Suppl. 
material 1: table $3. Boxplots indicate the mean values of each 
trait, quartiles and the range of the data (whiskers). Individual 
data points are overlaid as scatterplots. Letters indicate the signifi- 
cant differences among groups (PWRST with corrected p-values 
< 0.05). figure S3. Scatterplots of significant AIC environmental 
variables. For all female allopatric populations: ovipositor length 
vs. latitude (a.), elevation (b.), average soil percent clay (c.), mini- 
mum temperature (d.), and maximum temperature (e.). For all 
female allopatric populations: pronotum width vs. latitude (f.), 
longitude (g.), average soil % sand (h.), and minimum tempera- 
ture (i.). For female allopatric and sympatric Connecticut popula- 
tions: ovipositor length vs. minimum temperature (j.), maximum 
temperature (k.), and average soil % organic matter (I.). 
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